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The (l+l)-dimensional kinetic model of crystal growth with simulated self- attraction and random 
sequential or parallel dynamics is introduced and studied via Monte- Carlo simulations. To imitate 
the attraction of absorbing atoms the probability of deposition is chosen to depend on the number 
of the nearest-neighbor atoms surrounding the deposited atom so it increases with this number. As 
well the evaporation probabilities are chosed to roughly account for this self-attraction. The model 
exhibits the interface depinning transition with KPZ-type roughness behavior in the moving phase. 
The critical indices of the correlation lengths are = 0.82 ± 0.03, iy± = 0.55 ± 0.02 and the critical 
index of the growth velocity is 1.08 it 0.03 indicating the new universality class of the depinning 
transition. The critical properties of the model do not depend on the type of dynamics implemented. 

PACS numbers: 05.40.-a, 05.70.Fh 



I. INTRODUCTION 

The atomistic kinetic models of crystal growth studied extensively in recent years exhibit the interesting phenomenon 
of depinning transition between the phases of pinned and moving interface Q-Q. The interface roughness also changes 
drastically at this transition; it is finite in the pinned phase while it diverges in the moving phase of an infinite sample. 
The statistics and kinetics of the interface roughness in the moving (rough) phase are universally described b^ 
continuous Kardar-Parisi-Zhang (KPZ) equation [6] (or its limiting case - Edwards- Wilkinson (EW) equation |7]) as 
it captures the most relevant features of the process forming the interface: surface tension and lateral growth [8]. 

Yet KPZ equation can not fully describe the depinning transition as it does not predict the velocity of the interface 
growth and the roughening transition appears in it only in spatial dimensions d > 2 [1]. So the adequate description 
of crystal growth warrants the studies of the atomistic kinetic models of Sold-On-Solid (SOS) type devoid of these 
drawbacks. There are a number of such models imitating surface tension and lateral growth in various ways: via edge 
of plateau evaporation [3], polynuclear growth [4] or diffusion [5*]. Thus introduced correlations reduce the absolute 
value of interface roughness and make theoretical results qualitatively similar to the real picture of crystal growth. 

At the same time the various SOS models have rather different characteristics of the depinning transition such 
as critical indices of interface velocity and correlation lengths. Now it is not known how many universality classes 
of depinning transition can exist. To elucidate this point the studies of other possible SOS models are needed. So 
here we present one more model which can give a realistic picture of crystal growth. Its motivation relies on the 
fact that the major role in deposition-evaporation processes belongs to the self-attraction of the atoms of growing 
crystal. To account for this basic feature the deposition probability of the model should depend on the number of 
nearest-neighbor atoms surrounding the deposition cell being larger for more populated neighborhood. Apparently, 
the reverse situation should hold for the evaporation probabilities. 

Here we consider the Id model with such probabilities of deposition-evaporation processes and random sequential 
or parallel dynamics. We show that such simulation of self-attraction is sufficient to give a reasonably adequate 
description of depinning transition. The manifestation of this is, in particular, that the model has the roughness 
behavior of KPZ type. The critical behavior of the interface velocity and critical indices of the correlation lengths 
are shown to differ from that of the previously studied models which indicates the new universality class of depinning 
transition. 



II. MODEL 



The crystal growth from the non-crystalline (gaseous or liquid) phase is considered as the result of the deposition 
and evaporation processes. In the first one the randomly wandering crystal-forming atoms of the non-crystalline 
phase can occur near the crystal surface and stick to it, while in the second one the atoms evaporate from the crystal 
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FIG. 1: Three possible configuration of nearest sites around tfie deposition cell with three, two and one neighboring atoms. X 
denotes the empty cell and O shows the atom of crystalline phase. 

surface. The deposition cell in d = 1 case can have from one to three nearest neighbors as Fig. [T] shows. Accordingly 
to imitate the attraction of the deposited atom to those already belonging to the crystalline phase there should be 
three different probabilities of absorption Wk^ {k = 1^2,3 being the number of nearest neighbors in crystalline phase) 
such as 

1 > Ws > W2 > wi > 0. (1) 

Obviously, if k;i = 1 then ws = W2 = I and if ws = then wi = W2 = 0. Here we choose the simplest 
parametrization, complying with the condition ([1]) 

«;fe=/-^ 0<p<l (2) 

thus reducing the parameter space. 

To define the evaporation probabilities Uk in configurations of Fig. [T]we note that for /c = 3 the evaporating atom 
has always three nearest atoms, for /c = 2 it has two or three neighbors and for A: = 1 their number can be one, two 
or three. Then the assumed self-attraction can be roughly simulated by the condition 

1 > Ui > U2 > U3> 0. (3) 

The simple way to comply with it is to assume that no height changes occur with the equal probability q at all k. 
Then we have the evaporation probabilities 

Uk = l-q-Wk (4) 
for the configurations in Fig. [TJ The positivity of Uk imposes the limitation on the p, q values: 

p + g < 1. 

Thus the model assumes that the larger number of neighboring atoms the more effectively crystal retains the deposited 
atom among them and the evaporation takes place more often from the less populated (on the average) neighborhood. 
Just this one can expect as the result of the interatomic attraction. So the atoms are mainly absorbed and retained 
at the interface dips while their deposition on the local tops is more rare but the evaporation is more often. On the 
global scale this acts as the surface tension and induces the lateral growth of interface which result in smoothing 
its profile which withstand the randomness-induced roughness. Thus we may expect that the model captures the 
essential features of crystal growth and can give the adequate description of this process. 

Apparently there is a vast choice of probabilities complying with Eqs. ([H [3]) other than adopted here and some 
physical arguments can be used to define them. In particular, one may define the evaporation probabilities to depend 
strictly on the nearest-neighbor number as in the case of the absorption ones. Yet it needs the expansion of the update 
rules to involve the next-nearest neighbors of the upper atoms of crystal phase. This seems to be not quite necessary 
for our aim to demonstrate the principal validity of the crystal-growth mechanism based on the self-attraction of the 
crystal- forming atoms. 
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FIG. 2: The interface profiles and the average heights at ^ = 0.1 and p = 0.67 (a, b), p = 0.7 (c, d) at different times (1000 
-blue, 2500 - green, and 5000 - red MC steps). 

One can choose either random sequential or parallel updates for the model defined by Eqs. ^ [3j). Yet random 
sequential dynamics is more adequate as then at every time step we deal with definite interface profile while in the 
parallel updates the interference of changes in the nearby sites can not be generally taken into account. Still the 
alternative parallel updates of the odd and even sites can be implemented in which no interference appears as there 
are no correlations between the next-nearest sites in the update rules. This type of dynamics should give the results 
identical to those of random sequential one and our Monte-Carlo simulations corroborate this expectation. Moreover, 
the critical behavior of the model stays the same under the fully parallel updates; the slight changes appear only in 
the phase diagram. All the results presented below are obtained with the alternative parallel updates of the odd and 
even sites. 



III. MONTE-CARLO SIMULATIONS 



The Monte-Carlo simulations were performed on the samples of L = 100 — 400 sites with periodic boundary 
conditions for up to 10^ time steps. The interface profiles emerging at (7 = 0.1 and p = 0.67, 0.7 and at different times 
(1000, 2500, and 5000 MC steps) are shown in Fig|2l The appearance of the depinning transition is clearly seen in 
it. Quantitatively it is manifested in the appearance of the linear in time growth of the average interface height (cf. 

Fig. m,m) 
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FIG. 3: .Phase diagram. Circles are the results of MC simulations. 



Here hn{t) is the crystal height of site n at time t. The transition line p = Pc{q) is shown in Fig. [3l with good precision 
it can be approximated by the straight line 

Pc{q) = 0.715 - 0.387(7 

The two examples of the time dependence of the interface roughness 

L 

w{t)=L-'^\hn {t)-hav {t)\ 
n=l 

in the moving (rough) phase are shown in Fig. [H It shows that initial growth of w{t) complies with power low, 
w{t) . The average growth exponent [3 in it is found to be /3 = 0.35 ± 0.04 for q = 0,0.1 which is close to the 
KPZ value (3k pz = 1/3. To find the dynamical critical index z the simulations with the initial inclined interface 




FIG. 4: The time dependence of roughness for L = 200 in the moving phase (a) p — 0.7, = 0, (b) p = 0.7, q = 0.1. Straight 
lines correspond to the power law with the exponent 0.333. 

hn{t = 0) = n were performed in the moving phase. In this case initial w{0) = L/4 decays approximately exponentially 
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with time to the roughness' saturation value, w (oo) ~ <C (0) , a = /3z < 1, see Fig. O 

w{t)^w (0) exp {-t/r) + w (oo) (5) 
Here the usual scaling suggests [l| , Q that r ^ sufficiently close to the transition point where spatial correlation 
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FIG. 5: (a) Evolution of the interface with initial slope, q = 0.1, p = 0.8, L = 300. The profiles after t = 100, 3000, 6000 steps 
are shown; {h)w{t), solid line corresponds to Eq.©. 

length is much greater than L. These simulations at g = 0.1 andp = 0.8 for L = 100, 200, 300, 400 give z = 1.56±0.07 
(Fig. E^) which is close to KPZ value zrpz = 3/2 [l|, ^. For g = 0,p = 0.8 we get z = 1.96±0.1 (Fig. [Dd) indicating 
the crossover to the EW behavior with z = 2 ai p > 0.8. The r values in Fig. [6] are the average over 50-100 trials 
which appears to be sufficient for the saturation of the standard deviation. 




30 



W 20 



JO 



(b) 



























uia? i^io* 3x101* 
f 



FIG. 6: Double logarithmic plot r vs pc 
Straight line corresponds to z = 1.96. 



-p. (a) g = 0.1, p = 0.8. Straight line corresponds to z = 1.56, (b) g = 0, p = 0.8. 



We can not reliably conffim that at g = in the rough (moving) phase the EW behavior changes to the KPZ one 
nearer the transition point as it needs the studies of much larger samples. We do not also observe the region of the 
EW behavior w {t) ~ t^^^ at t less than some tc 0, cf. Fig. [4^. This is probably due to a small tc value. Yet this 
crossover can be detected in the critical index of the longitudinal (time) correlation length at different p in the 
pinned phase. To find it we again determine the characteristic time r for the vanishing of initial slope but now as a 
function of p — Pc Fig. [71 illustrates the results of such simulations. The expected scaling behavior r ~ {pc — p)^" 
holds for ^ = 0.1 at ah 0.5 < p < = 0.68 with v\\ = 0.82 ± 0.03, cf. Fig. But for q = the crossover from the 
mean-field (EW) i^y = 1 to the u\\ = 0.83 ib 0.04 takes place at p = 0.66, see Fig. [8]3. So assuming that in the critical 
regions of both phases z/y and z are the same we may suppose that at = the critical behavior of the model also 
belongs to that of KPZ class. Hence the index of the spatial correlation length is u± = uu/z = = 0.55 ± 0.02. 
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FIG. 7: (a) Evolution of the interface with initial slope, g = 0, p = 0.68, L = 100. The profiles after t — 150, 750, 1500 steps 
are shown; (b) w{t)^ intersection of straight line with the time axis defines the characteristic decay time. 
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FIG. 8: (a)Double logarithmic plot r vs L. (a) ^ = 0.1, L = 100. Straight line corresponds to z/y 
Straight Unes correspond to = 0.83 and =1. r is the average over 50-100 trials. 



0.82; (b) ^ = 0, L = 100. 



It is just slightly above the mean- field value = z = 1/2. These u± and u\\ values are rather low as compared to 
those known for other (1+1) models. Thus in models of Refs. Q they are equal to those of directed percolation 
i^±,DP ^ 1.1 and J^w^DP ^ 1-73. 

For q = and q = 0.1 we have also determined in the moving phase the growth velocity index 7 in the relation 
V ^ {p — PcV . In both cases it appears to be close to unity, 7 = 1.08 ± 0.03, see Fig. [9l It differs significantly from 
that of the models in which it can be identified with the index of the longitudinal correlation length of the directed 
percolation, z^y 2:)p 1.73 0, [1]. 

IV. CONCLUSIONS 



The present results show that simulation of the attraction between the crystal-forming atoms in the kinetic growth 
model can provide the realistic picture of crystal growth even in the absence other mechanisms withstanding the 
interface roughening such as diffusion. The roughening transition in the model belongs to the KPZ class which can 
be changed to EW one farther from transition point. The depinning transition in it differs in the critical behavior 
from the previously studied models by the lower indices z^y , and 7. 

We should also point out that the model lacks the tilting symmetry which KPZ equation has [8], the initial slope 
in it always vanishes and the interface becomes parallel (on the average) to the initial surface, see Figs. O H It 
seems probable that other SOS models [l]-[5] also have the preferable growth direction. Meanwhile this property 
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FIG. 9: Double logarithmic plot v vs p — pc- (a) g = 0, L = 200. Straight line corresponds to 7 = 1.08; {h)q = 0.1, L = 200. 
Straight line corresponds to 7 = 1.08. 

becomes lost in the renormalization group procedure for the general macroscopic growth equations with (Gaussian) 
noise leading invariably to the KPZ equation on the large space and time scales [l|, [H. This may point out that 
the Langevin-type stochastic equations can not adequately describe the variety of the depinning transitions and some 
other more general approaches to their macroscopic characterization should be developed. 
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